function [equalorder, ...
          Ugridreduced, sgreduced,...
          id_Ugridreduced,...
          s_id_gridreduced, sU_id_gridreduced,...
          n_U_sets, Tcoord, indices_pairs,...
          numeqconstraints,coordrowseq, fillAeq, beq]=useful_anyparam1_ind(u01grid, u11grid, u21grid,u02grid, u12grid, u22grid, sg,PYX_cond)     
%% U sets
n_U_sets=3; %number of U-sets (number of dimensions of the distribution of the taste shock differences)
Ugrid=cell(1,n_U_sets); %cell 1xn_U_sets
Ugrid{1}=[u01grid-u11grid u02grid-u12grid Inf*ones(sg,1)  -Inf*ones(sg,1)]; 
Ugrid{2}=[u01grid-u21grid u02grid-u22grid Inf*ones(sg,1)  -Inf*ones(sg,1)]; 
Ugrid{3}=[Inf*ones(sg,1) -Inf*ones(sg,1)  u21grid-u11grid  u22grid-u12grid]; 

%% Reduce number of grid points to check by applying our partitioning arguments
%equalorder (sgx1) assigns the same index to parameters vectors such that Ugrid{1}, Ugrid{2}, ..., Ugrid{n_U_sets} have the same order

%%%%%%%% Criterion 1 %%%%%%%%
d=cell(n_U_sets,1);
index=cell(n_U_sets,1);
for j=1:n_U_sets 
    [sortedUjgrid,index{j}] = sort(Ugrid{j},2); 
    %sortedUjgrid gives Ugrid{j} with each row sorted from smallest to largest; it is sgxn_columns_Ujgrid
    %for each row i of sortedUjgrid and for k=1,2,...,n_columns_Ujgrid, index{j}(i,k) gives the column position in Ugrid{j}(i,:) of sortedUjgrid(i,k); it is sgxn_columns_Ujgrid
    d{j} = diff(sortedUjgrid,[],2) == 0; %sg x (n_columns_Ujgrid-1)
    %Suppose that Ugrid{j} has 3 columns. Then, for each row i of sortedUjgrid, d{j} will be: 
    %[1 1] if sortedUgridj(i,2)-sortedUgridj(i,1)=0 and sortedUgridj(i,3)-sortedUgridj(i,2)=0
    %[1 0] if sortedUgridj(i,2)-sortedUgridj(i,1)=0 and sortedUgridj(i,3)-sortedUgridj(i,2)~=0
    %[0 1] if sortedUgridj(i,2)-sortedUgridj(i,1)~=0 and sortedUgridj(i,3)-sortedUgridj(i,2)=0
    %[0 0] if sortedUgridj(i,2)-sortedUgridj(i,1)~=0 and sortedUgridj(i,3)-sortedUgridj(i,2)~=0
end
%%%%%%%% Criterion 2 %%%%%%%%
n_columns_Ujgrid=4;
Ugrid_sum=zeros(sg,n_columns_Ujgrid+n_columns_Ujgrid^2);
for g=1:sg
Ugrid_sum(g,:)=[Ugrid{1}(g,:) reshape((Ugrid{2}(g,:)+Ugrid{3}(g,:).')',[],1)'];
end

[sortedUgrid_sum,index_sum] = sort(Ugrid_sum,2); 
d_sum = diff(sortedUgrid_sum,[],2) == 0; %sg x (n_columns_Ujgrid*n_U_sets-1)

%%%%%%%% Combine %%%%%%%%
[~,~,equalorder] = unique([index_sum d_sum [index{:}] [d{:}] ],'rows', 'stable'); %sgx1 

%Reduced grid of parameter vectors and U sets
[~,index_final,~]=unique(equalorder, 'stable'); %one row from equalorder per ordering group together with its index
Ugridreduced=cell(1,n_U_sets);
for j=1:n_U_sets 
    Ugridreduced{j}=Ugrid{j}(index_final,:); %sgreduced x n_columns_Ujgrid
end
sgreduced=size(Ugridreduced{j},1); %size of the reduced grid of parameter values


%% Construct indices to keep track of equal elements inside Ugridreduced{j} for j=1,...,n_U_sets
id_Ugridreduced=cell(n_U_sets,1);
for j=1:n_U_sets 
    id_Ugridreduced{j} = NaN(size(Ugridreduced{j})); % %sgreduced x n_columns_Ujgrid
    for k = 1:sgreduced
        [~, ~, id_Ugridreduced{j}(k,:)] = unique(Ugridreduced{j}(k,:), 'stable'); 
    end
end

s_id_gridreduced=zeros(sgreduced, n_U_sets); %useful to set number of unknowns of the linear programming
for j=1:n_U_sets
    s_id_gridreduced(:,j)=max(id_Ugridreduced{j},[],2); %number of distinct elements of Ugridreduced{j}
end
sU_id_gridreduced=prod(s_id_gridreduced,2); %sgreducedx1 reporting the cardinality of the Cartesian product of the U sets
                                            %sU_id_gridreduced(g) is the number of unknowns of the LP for each g-th parameter value
s_id_gridreduced_2=ones(1, n_U_sets)*n_columns_Ujgrid; %number of elements of Ugridreduced{j} (count as separate even equal elements)
sU_id_gridreduced_2=prod(s_id_gridreduced_2); %cardinality of the Cartesian product of the U sets     

%% Invariant components of the linear programming problem 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%% Inequality constraints %%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

vectors = arrayfun(@(x) {1:x}, s_id_gridreduced_2); 
n = numel(vectors);
Tcoord_temp = cell(1,n);
[Tcoord_temp{:}] = ndgrid(vectors{:}); 
Tcoord_temp = cat(n+1, Tcoord_temp{:});
Tcoord = reshape(Tcoord_temp,[],n); %matrix of dimension sU_id_gridreduced_2xn_U_sets reporting the column coordinates of all the possible n_U_sets-tuples from the U-sets
%For example if n_U_sets=3 and s_id_gridreduced_2 is [2 3 1], then 
%[ca, cb, cc] = ndgrid(1:2, 1:3, 1:1);
%Tcoord = [ca(:), cb(:), cc(:)];   
%Tcoord= [1     1     1
%         2     1     1
%         1     2     1
%         2     2     1
%         1     3     1
%         2     3     1];
indices_pairs=pairIndices(sU_id_gridreduced_2); %row indices of all possible unordered pairs of rows from the matrix obtained by taking the Cartesian product of the U sets
                                         %[sU_id_gridreduced_2*(sU_id_gridreduced_2-1)/2]x[2]   
%Continuing the example above
%indices_pairs=[1 1; 1 2; ...; 1 27; 2 3; ...; 2 27; ... 26 27]

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%% Equality constraints %%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

numzeros=4*3+3*3+4*4;
numeqconstraints=6+(numzeros)+1;

coordrowseq=[1 1 1 ... %observational equivalence
             2 2   ...
             3     ...
             4 4 4 ...
             5 5   ...
             6     ...
             7:1:6+numzeros ... %zeros
             6+numzeros+1]; %ones;  
                         
fillAeq=[1 -1 -1 ... %observational equivalence
         1 -1    ...
         1       ...
         1 -1 -1 ... %observational equivalence
         1 -1    ...
         1       ...
         ones(1,numzeros) ... %zeros
         1]; %ones

beq=[ PYX_cond(1)-1;... %observational equivalence
      PYX_cond(2);...
      PYX_cond(3);...
      PYX_cond(4)-1;... %observational equivalence
      PYX_cond(5);...
      PYX_cond(6);...
      zeros(numzeros,1);... %zeros
      1];%ones  %[#equalityconstraints]x1      
end